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We develop the theory of canonical-dissipative systems, based on the assumption that both 
the conservative and the dissipative elements of the dynamics are determined by invariants of 
motion. In this case, known solutions for conservative systems can be used for an extension 
of the dynamics, which also includes elements such as the take- up/dissipation of energy. This 
way, a rather complex dynamics can be mapped to an analytically tractable model, while still 
covering important features of non-equilibrium systems. 

In our paper, this approach is used to derive a rather general swarm model that considers 
(a) the energetic conditions of swarming, i.e. for active motion, (b) interactions between the 
particles based on global couplings. We derive analytical expressions for the non-equilibrium 
velocity distribution and the mean squared displacement of the swarm. Further, we investigate 
the influence of different global couplings on the overall behavior of the swarm by means of 
particle-based computer simulations and compare them with the analytical estimations. 
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I. INTRODUCTION 

The collective motion of biological entities, like schools 
of fish, flocks of birds, herds of hoof animals, or swarms 
of insects, recently also attracted the interest of physi- 
cists. Here, the question how a long-range order between 
the moving entities can be established is of particular in- 
terest. Consequently, some of the more biologically cen- 
tered questions of swarming behavior, namely about the 
resons of swarming or the group size dependence, have 
been dropped so far in physical swarm models. The main 
focus was rather on the emergence of coherent motion in 
a "swarm" of locally or globally coupled particles. 

For the coupling different assumptions have been pro- 
posed (cf. also Sect. IV), such as the coupling of the par- 
ticles' individual orientation (i.e. direction of motion) to 
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the mean orientation of the swarm 0, |^ , or the coupling 
of the particles' individual position to the mean position 
(center of mass) of the swarm On the other hand, also 
some more local couplings have been considered, such as 
the coupling of the particle's individual velocity to a local 
average velocity ^ |5[ ||]. A different class of models 
further assumes a local coupling of the particles via a self- 
consistent field that has been generated by them |Q, |j. 
This models the case of chemical communication between 
the particles as widely found in biology. For example, the 
streaming behaviour of myxobacteria |l^, |ll| , or the 
directed motion of ants |lj, ^ Q or the coherent move- 
ment of cells have been described and modeled by 
means this kind of local coupling. 

But long-range or short-range coupling of the particles 
is only one of the prerequisites that account for swarm- 
ing. Another one is the active motion of the particles. Of 
course, particles can also move passively, driven by ther- 
mal noise, by convection, currents or by external fields. 
This kind of driving force however does not allow the 
particle to change its direction of motion, or velocity 
etc. itself. Recent models of self-driven particles which 
are used to simulate swarming behavior ||^, |l6[ |l^ usu- 
ally just postulate that the entities move with a certain 
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non-zero velocity, without considering the energetic im- 
plications of active motion. In order to do so, we need 
to consider that the many-particle system is basically an 
open system which is driven into non-equilibrium. 

To this end, our approach to swarming is based on 
the theory of canonical- dissipative systems. This the- 
ory - which is not so well-known even among experts 
- results from an extension of the statistical physics 
of Hamiltonian systems to a special type of dissipa- 
tive systems, the so-called canonical-dissipative systems 
0, ^ 1^, H, H H. The term dissipative means 
here that the system is non-conservative and the term 
canonical means, that the dissipative as well as the con- 
servative parts of the dynamics are both determined by 
a Hamilton function H (or a larger set of invariants of 
motion, see Sect. || for details). 

This special assumption allows in many cases exact 
solutions for the distribution functions of many-particle 
systems, even in far-from-equilibrium situations. This 
was known already to pioneers as Poincare, Andronov 
and Bautin who gave exact solutions for a special class 
of nonlinear oscillators In recent work the 

properties of canonical-dissipative systems were used to 
find exact solutions for dissipative solitons in Toda lat- 
tices which are pumped with free energy from external 
sources. The starting point was the well-known Toda the- 
ory of soliton solutions in one-dimensional lattices with 
a special nonlinear potential (exponential for compres- 
sion and linear for expansion) of the springs p7| , p8| , 
which was then extended towards a non-conservative sys- 
tem ^ |2^. In another work an application to 
systems of Fermi and Bose particles was given. 

In this paper, we will apply the theory of canonical- 
dissipative systems to the dynamics of swarms. We start 
with an outline of the general theory in Sect. ||, describ- 
ing first the deterministic and then the stochastic ap- 
proach. As a first step towards the dynamics of swarms, 
in Sect. ^ we investigate the energetic conditions of 
swarming, i.e. we discuss the conditions for active mo- 
tion and their impact on the distribution function and 
the mean squared displacement of the particles. In Sect. 
IV we introduce global interactions between the particles 
which may account for swarming, i.e. for the maintenance 
of coherent motion of the ensemble. The global inter- 
actions are choosen with respect to the general theory 
of canonical-dissipative systems. By means of computer 
simulations we show how different kinds of coupling af- 
fect the overall behavior of the swarm. In Sect. ^ we 
conclude with some ideas how to generalize the approach 
presented in this paper. 



II. GENERAL THEORY OF 
CANONICAL-DISSIPATIVE SYSTEMS 



A Dynamics of Canonical-Dissipative Systems 

Let us consider a mechanical many-particle sytem with 
/ degrees of freedom i = 1, / and with the Hamilto- 
nian H{qi ... 9/, pi ■ ■ ■ Pf) The corresponding equations 
of motion are 



dqi dH dpi dH 
dt dpi ' dt dqi 



(1) 



Each solution of the system of equations (^ 

Pi ^Pi{t) ; qi = Qiit) (2) 

defines a trajectory on the plane H — E =const. This 
trajectory is determined by the initial conditions, and 
also the energy E = H{t = 0) is fixed due to the initial 
conditions. We construct now a canonical-dissipative sys- 
tem with the same Hamiltonian, where g{H) denotes the 
dissipation function: 



dt dq, 'dp. 



(3) 



In order to elucidate this kind of canonical-dissipative 
dynamics, we will consider different examples for the 
dissipation function in the following. In general we will 
only assume that g{H) is a nondecreasing function of the 
Hamiltonian. The canonical-dissipative system, eq. (|^), 
jTs], 11^, |23| does not conserve the energy because of the 
following relation: 



-dt=-'^''^2:[9f^ 



(4) 



Whether the total energy increases or decreases, conse- 
quently depends on the form of the dissipation function 
g{H). In the simplest case, we may consider a constant 
friction: 



g{H) = 70 > 



(5) 



As far as g{H) is always positive, the energy always de- 
cays. 

In a more interesting case the dissipative function g{H) 
has a root for a given energy Ei: g{Ei) — 0. Let us fur- 
ther assume that at least in certain neighbourhood of Ei 
the function g{H) is increasing. With these assumptions 
the states with H < Ei are pumped with energy due to 
the negative dissipation, while energy is extracted from 
the states with H > Ei. That means, any given initial 
state with H{t = 0) < Ei will increase its energy until it 
reaches the shell H{t) = Ei, while any given initial state 
with H{t — 0) > El will decrease its energy until the 
shell H(t) = El is reached, too. Therefore the solution 
of eq. (|j) converges to the energy surface H = Ei. On 
this surface the solution of eq. (^ agrees with one of the 
possible solutions of the original Hamiltonian eq. (|l|) for 
H = Ei. 
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The simplest ansatz for g{H) with a root g{H) = is 
a linear dissipation function: 

giH) = C{H - El) (6) 

With respect to eq. and the discussion above, it is 
obvious for this case that the process comes to rest when 
the shell H{t) = Ei is reached. The relaxation time is 
proportional to the constant C~^. 

We note that the linear dissipation function eq. (||) 
has found applications in Toda chains Here, the fact 
that on the shell H = Ei the trajectory should obey the 
original conservative canonical dynamics, has been used 
to derive exact solutions for canonical- dissipative Toda 
systems. 

A rather general nonlinear and non-decreasing dissi- 
pation function, which has been proposed in |29 , reads: 



g{H) =70-71 



(1 + ^) 



l + Aexp{(3H) 



(7) 



Here 70 > represents the normal positive friction, 
7i > represents a kind of negative friction. A is a 
dimensionless constant and /3 is a parameter with the 
meaning of a reciprocal temperature. For 71 < 70 the 
friction is always positive, i.e. energy is extracted. For 
the opposite case, 71 > 70, we have negative friction and 
the system is pumped with energy at least in some parts 
of the phase space. This allows to drive the system into 
situations far from equilibrium, therefore in the following 
we may assume 71 > 70 > 0. 

In the limit /3 — > and A ^ the dissipative function 
eq. (0) reduces to the linear case, eq. (||), discussed above. 
On the other hand, for small values of /3 and finite A we 
get 



g{H) =70-71 



(1 + ^) 



1 + A + Af3H 



(8) 



which yields g{H) (70 — 71) < for H ^ 0, and 
g{H) ^ 7o, if -ff ^ cxD. This way, the existence of a root 
g{Ei) = is always guaranteed. 

Before investigating a special case for g{H) in Sect. 
HI, let us introduce a generalization of the formalism. 
Instead of a driving function g{H) which depends only 
on the Hamiltonian we may include the dependence on 
a larger set of invariants of motion, Iq, Ii, I2, Is, for 
example: 

• Hamilton function of the many-particle system: 
Io = H 

• total momentum of the many-particle system: Ii = 



• total angular momentum of the many-particle sys- 
tem: I2 = L 

The dependence on this larger sets of invariants 
may be considered by defining a dissipative potential 



motion, eq. 



The canonical-dissipation equation of 
I) is then generalized towards: 



df^ ^ _dH_ _ dG{Io,h,l2,-) 
dt dqi dpi 



(9) 



Using this generalized canonical-dissipative formalism, 
by an appropriate choice of the dissipative potential G 
the system may be driven to particular subspaces of the 
energy surface, e.g. the total momentum or the angu- 
lar momentum may be prescribed. Different examples for 
this will be also discussed in Sect. IV. 
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Stochastic Theory of Canonical-Dissipative 
Systems 



We will now investigate an approach to the station- 
ary probabilities which is based on Langevin and Fokker- 
Planck equations. The Langevin equations are obtained 
by adding a white noise term Ciit) to the deterministic 
eq. (|): 

The essential assumption is that both the strength of the 
noise, expressed in terms of D{H) and the dissipation, 
expressed in terms of the friction function g{H) depend 
only on the Hamiltonian H. 

With respect to eq. ([l0|), the corresponding Fokker- 
Planck equation for the probability distribution 
p(qi...qf ,pi...pf) of the many-particle system reads: 



dp \ - dp 



dqi 



^ dH dp 
^ dpi dpi 



E 



d_ 

dpi 



g{H)p + D{H) 



dp 

dpi 



An exact stationary solution of eq. (|llD reads: 



H 



p°{qi...qf,pi...pf) ^ Q ^ exp J dH 



D{H') 



(11) 



(12) 



The derivative of p° vanishes if g{H) = 0, which means 
that the probability distribution is maximal at the sur- 
face H = El. 

For the special case of a linear dissipation function, eq. 
(^), we find the stationary solution: 

p0(gi...g/,pi...p/) = exp ( ''^ ^^fj^~^^ ) (1^) 

In the limit of very strong pumping, iJ 3> i?i , this prob- 
ability distribution reduces to a kind of microcanonical 
ensemble corresponding to the energy Ei . 
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The existence of exact solutions for the probabiHty 
distribution, eq. (^2|), allows to derive several thermo- 
dynamic functions for the many-particle system, such as 
the mean energy: 



H 



U 



Q-^ j dH J{H) Hcxp J dH 



Dm 



where the Jacobian J{H) is defined by: 

dqi...dqf dpi...dpf — J{H) dH 



(14) 



(15) 



The entropy follows from the Gibbs formula which yields 
here: 

S = +kB \nQ-kB J dH J{H) x 



H 



Q ^ exp 



9{H) 
D{H) 



H 



dH 



Dm 



(16) 



Further, the system has a Lyapunov functional K which 
is provided by the Kullback entropy: 



^ J dqi...dqfdpi...dpf pin 
This gives explicitely: 



(17) 



K[p,po] = InQ - — 
Kb 




(18) 



This functional is always nonincreasing. 

If we want to generalize the description by means of 
the dissipative potential G{Io, Ii, h, ■■■) as used in eq. 
(^, the situation is more difficult. But at least if the 
noise strength D{H) is a constant {D =const.), the cor- 
responding Fokker-Planck equation 



dp 

dt 
d 



Pi 



dp 



^ dH dp 
^ dpi dpi 



^ dp^ 



dG{h,h,l2,-) 
dp., ' 



D 



dp_ 

dpi 



(19) 



is still solvable. Due to the invariant character of the 
/fe, the l.h.s. of eq. ( |l9| ) disappears for all functions 
of Ik- Therefore we have to search only for a function 
G{Io, Ii, I2, ■■■), for which also the collision term disap- 
pears. This way we find the stationary solution of eq. 



p^{qi...qf,pi...pf) = 



exp 



(20) 

The derivative of po vanishes if G{Iq, Ii, I2, ■■) — min. 
which means that the probability is maximal on the at- 
tractor of the dissipative motion. 



III. ENERGETIC CONDITIONS OF 
SWARMING 

In this section, we want to apply the general descrip- 
tion outlined above to swarm dynamics. Basically, at 
a certain level of abstraction a swarm can viewed as 
a many-particle system with some additonal coupling 
which would account for the typical correlated motion 
of the entities. In addition, also some energetic condi- 
tions must be satisfied in order to keep the swarm mov- 
ing. Thus, active motion of the particles with a non-zero 
velocity is another basic ingredient for swarming. Since 
the conditions for active motion are dropped in many of 
the currently discussed swarm models ^, |l^, we first 
want to investigate the energetic conditions of swarm- 
ing, before turning to the second ingredient, coupling of 
individual motion. 



A Conditions for Active Motion 

If we omit interactions between the particles, the 
Hamiltonian of the many-particle system is of the simple 
form: 



N N 2 

H^yH, = yp^ 



(21) 



This allows to reduce the description level from the A^- 
particle distribution function to the one-particle distri- 
bution function: 



p{qi...qf,pi...pf) = J|p(ri,pi) 



(22) 



where the variable r j is now used for the space coordinate 
of the particle z, and Pj stands for the momentum or, if 
m = 1 is used in the following, for the velocity of the 
particle, respectively. The center of mass of the swarm is 
defined as 



(23) 



whereas the mean momentum P and the mean angular 
momentum L are defined as: 

P-l^T.P^ (24) 

^=^E^» = ^E^'XP^ (25) 

With g{H) ~ g{Pi) and D = const., the Langevin eq. 
( |To| ) for each particle i reads in the absence of an external 
potential: 



= Pi ; P^ 



-g{pf)p, + {2DY^'m (26) 



where the strength D of the stochastic force results from 
the Einstein relation: 



(27) 



5 



T is the temperature and fc^ is the Boltzmann constant. 

For the dissipation function gijP'), we use again the 
general nonhnear ansatz of eq. (Q) in the hmit of smaU 
values of [3 and finite A, eq. (jg). By means of the substi- 
tutions 



7i(l + ^) =d2So ; f3A = 2d2 (28) 



(1 + ^) 



and with H of eq. (|2l|), eq. (H) can be written in the form: 



9{P^) = 7o 



d2P^ 



(29) 



We note that all noninteracting systems with g = g^p^) 
are of canonical-dissipative type. 

Eq. ( p9[ ) agrees with the velocity dependent non-linear 
friction function previously used in a model of active 
Brownian particles These are driven Brownian 

particles, which move due to the influence of a stochastic 
force, but additionally are pumped with energy due to 
a velocity-dependent dissipation function g{p^), eq. ( p9| ) 
which is plotted in Fig. 0. 
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FIG. 1: Velocity-dependent dissipation function, g{p^), 
eq. (29) vs. p. The velocity ranges for "pumping" {g{p^) < 



0) and "dissipation" {g{p^) > 0) are indicated. Parame- 
ters: So = 10.0, 70 = 20.0, c = 1.0, d2 = 10.0. 

The second term of the r.h.s. in eq. (p9| ) which results 
from the pumping of energy, has been physically sub- 
stantiated in our earlier work We have assumed 
that the Brownian particles are able to take up energy 
from the environment at a constant rate so: which can 
be stored in an internal depot, e. The internal energy can 
be converted into kinetic energy with a velocity depen- 
dent rate d{p) — d2P^, which results in an additional 
acceleration of the Brownian particle in the direction of 
movement. The value of the internal energy depot may 
be further decreased due to internal dissipation processes, 
described by the constant c. Thus the resulting balance 
equation for the energy depot reads: 



de 

— = So - c e{t) - d2 p^ e{t) 



(30) 



If we assume that the internal energy depot relaxes fast 
compared to the motion of the particle, we find the qua- 
sistationary value: 



eo 



So 



C -I- C?2P^ 



(31) 



which eventually leads to the second term of the r.h.s. in 



eq. 

Dependent on the parameters 70, ^2, soi c, the dissipa- 
tion function, eq. (p9[), may have a zero, where the friction 
is just compensated by the energy supply. It reads in the 
considered case: 



£0 
70 



c 

d'2 



(32) 



We see that for p < pq, i.e. in the range of small mo- 
mentums pumping due to negative friction occurs, as an 
additional source of energy for the Brownian particle. 
Hence, slow particles are accelerated, while the motion 
of fast particles is damped. 

For sod2 < 7oc, we find no real-valued root of eq. (|3^). 
This is the case of subcritical pumping, where the particle 
will move more or less like a simple Brownian particle. 
However, provided the existence of a non-zero momentum 
Po, i.e. for a supercritical pumping, the particle will be 
able to move in a "high velocity" or active mode H) , 
which displays several non-trivial features of motion, as 
will be shown by means of computer simulations in the 
next section. 



B Distribution Function and Dissipative Potential 

Due to the pumping mechanism discussed above, the 
conservation of energy clearly does not hold for the 
particle, i.e. we now have a non-equilibrium, canonical- 
dissipative system as discussed in Sect. |^ This results in 
deviations from the known Maxwellian velocity distribu- 
tion of an equilibrium canonical system. 

As pointed out in Sect. |IIB| , the probability density 
for the velocity p{p,t) obeys the Fokker-Planck eq. (|Tl|), 
which reads for the special case of the dissipation func- 
tion, eq. rt29) , and in the absence of an external potential: 



dt 



d_ 
dp 



70 



d2So 
+ d2P^ 



pp{p,t) + D 



dp{p,t) 



dp 
.(33) 

We mention that Fokker-Planck equations with nonlinear 
friction functions are discussed in detail in | |3^ . 

The stationary solution of eq. (|3^ ) is given by eq 
which reads in the considered case explicitely: 



P°(P) 



Co exp 
Co I 1 



D 

d2P^ ^ ^ 



exp 



70 2 
'2dP 



(34) 
(35) 
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where Cg results from the normahzation condition. 
Go(p^) is the special form of the dissipation potential 
G(/o,/i,/2, ■■•) considering only Iq = H as the invariant 
of motion, and further H as given by eq. (|2|). It reads 
explicitcly: 

„2 / ^^p2^ 



Gilo) = Goip' 



p 

7oy 



So , 



1 



(36) 



Compared to the Maxwellian velocity distribution of 
"simple" Brownian particles, a new prefactor appears 
now in eq. ( ^ ) which results from the additional pump- 
ing of energy. For a subcritical pumping, 591^2 < c7o, 
where we do not find a real-valued root of the dissipation 
function, eq. (|29|), only an unimodal velocity distribution 
results, centered around the maximum pp = 0. However, 
for supercritical pumping, 50^2 > C70, if the root of 7(p^) 
is real, we find a crater-like velocity distribution, which 
indicates strong deviations from the Maxwell distribution 

. . . . n . 

This is also shown in Fig. H which presents computer 

simulations of the velocity distribution of 10.000 particles 
after a sufficiently long time (only the x-dimension of the 
2-d simulation is shown). For the supercritical case, two 
distinct peaks of the velocity distribution are found at 
= {—0.63, +0.63}. The values of these maxima agree 
with the deterministic result for the stationary velocity, 
eq. (H). 

We note that non-Maxwellian velocity distributions for 
active motion have been also observed experimentally in 
cells, such as granulocytes p6|, p7|. 



C Mean Squared Displacement and Stationary 
Values 

As Fig. 1^ shows, the momentum distribution is cen- 
tered around p = both for subcritical and supercritical 
pumping. If we consider a nearly spherical swarm of par- 
ticles in the two-dimensional space r = {x, y} as in the 
computer simulations in this section, its center of mass, 
eq. (2^) and mean momentum, eq. ( [24| ) will come to rest. 
Thus they are not affected by the pumping, but other 
quantities are, such as the mean squared displacement: 

2\ 



AR^t) 




(37) 



In the limit of pure Brownian motion, it is known that 
the mean squared displacement increases in time as: 



AR^it) = ADrt 



(38) 



where Dr = ksT/^o = ^/7o is spatial diffusion co- 



efficient. Thus, eq. (^8|) will be the lower limit for sub- 
critical pumping of the particles. Contrary, in the case 
of supercritical pumping it has been shown [BSl p8| that 




FIG. 2: Velocity distribution p{px) for a swarm of 10.000 
particles aX t — 1.000 (i.e. in the stationary regime). 
Only of the 2-d simulation is shown, (top) subcriti- 
cal pumping ^2 — 1-0, (bottom) supercritical pumping, 
d2 — 10.0. Other parameters: D — lO^"^, sq = 10.0, 
7o = 20.0, c = 1.0. Initial conditions: r,(0) = {0.0,0.0}, 
p,(0) = {0.0, 0.0} for aU particles. 



the mean squared displacement will grow in time approx- 
imately as 



AR^{t) 



H 
D 



(39) 



where vq is given by eq. (p2|). Consequently, the diffusion 
coefficient in eq. ( |38| ) has for the case of supercritical 
pumping to be replaced by an effective spatial diffusion 
coefficient: 

2 



D 



off 



2D 



1 

2D 



£0 

70 



c 



(40) 



This result holds for noninteracting particles in the limit 
of relatively weak noise intensity D and/or strong pump- 
ing and will therefore give an upper limit for AR^ (t) . We 
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note the high sensitivity with respect to noise expressed 
in the scahng with (1/-D). 

Fig. ^ shows the mean squared displacement of a 
swarm of 2.000 particles both for the case of sub- and su- 
percritic al p umping together with the theoretical results 
of eqs. (p8[), (39). We see that for long times the com- 
puter simulations for supercritical pumping agree very 
well with eq. 




FIG. 3: Mean squared displacement Ai?2(i), eq. (|^ of 
a swarm of 2.000 particles as a function of time, (a) su- 
percritical pumping, ^2 — 10.0, (b) subcritical pumping, 
d2 — 1.0. D — 10~^, for the other parameters and the 
initial conditions see Fig. ^j. The additional curves give 
the theoretical results of eq. ( |39| ) (c) (upper limit) and 
eq. (Il) (d) (lower limit). 



Another quantity affected by the sub/supercritical 
pumping is the stationary velocity pg = Wq, eq. 
(^). In 2d, these stationary velocities define a cylin- 
der, Px + pI = Pq, in the four-dimensional state space 



{x, y,Px,Py} which attracts all deterministic trajectories 
of the dynamic system Fig. ^ shows the results of 
computer simulations iorp^{t) for the case of supercriti- 
cal pumping. The convergence toward the theoretical re- 
sult, eq. (p2|) can be clearly observed. 



IV. GLOBALLY COUPLED SWARMS 

So far we have neglected any coupling within the many- 
particle ensemble. This leads to the effect that the swarm 
eventually disperses in the course of time, whereas a 
"real" swarm would maintain its coherent motion. A 
common way to introduce correlations between the mov- 
ing particles in physical swarm models is the coupling to 
a mean value. For example, in ||l|, ^ the coupling of the 
particles' individual orientation [i.e. direction of motion) 
to the mean orientation of the swarm is discussed. Other 
versions assume the coupling of the particles' velocity to 
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FIG. 4: Averaged squared velocity p'^{t) — ^/NY^p1{t), 
eq. (p6|), of a swarm of 2.000 particles as a function of 
time. D = 10""*, for the other parameters and the ini- 
tial conditions see Fig. p|. T he dash-dotted line gives the 
stationary velocity, eq. (ft 21). 



a local average velocity^ which is calculated over a space 
interval around the particle ^ . 



A Coupling to the Center of Mass 

In this paper, we are mainly interested in global cou- 
plings of the swarm, which fit it into the theory of 
canonical-dissipative systems outlined in Sect. ||. As the 
most simple case we may first discuss the global coupling 
of the swarm to the center of mass, eq. (p3|). That means 
the particle's position is related to the mean position 
of the swarm R via a potential U{ri,R). For simplicity, 
we may assume a parabolic potential: 



U{r,,R)^^{r-Rf 



(41) 



The harmonic potential generates a force directed to the 
center of mass which can be used to control the dispersion 
of the swarm. It reads in the considered case: 



N 



VC/(r)=a(r,-i2) = -^(r, 

i=i 

With eq. (|2 
many-partic 



) , the corresponding Langevin eq. (|2 
e system reads explicitely: 



(42) 



of the 



N 



i2D)'/^m (43) 



Hence, in addition to the dissipation function there is 
now an attractive force between each two particles i and 
j which depends linearly on the distance between them. 



8 



With respect to the harmonic interaction potential eq. 
( ^ ) , we call such a swarm a harmonic swarm . 

Strictly speaking, the dynamical system of eq. ( [43|) is 
not a canonical-dissipative one, but as shown in it 
may be reduced to this type by some approximations, 
which will be also discussed below. We note that this 
kind of swarm model has been previously investigated in 
1^ for the one-dimensional case, however with a differ- 
ent dissipation function ^(p^), for w hich w e use eq. ( |29| ) 
again. Obviously, as shown in Sect. IIIB swarming will 
occur only for supercritical conditions. 

With the assumed coupling to the center of mass, R, 
the motion of the swarm can be considered as a superpo- 
sition of two motions: (i) the motion of the center of mass 
itself, and (ii) the motion of the particles relative to the 
center of mass. Taking into account that the noise acting 
on the different particles is not correlated, the center of 
mass for the assumed coupling obeys a force-free motion. 



1 ^ 

R^P, P^-j;^ T.9ip')P 



(44) 



Because of the nonlinearities in the dissipation func- 
tion g{p^) both motions (i) and (ii) cannot be sim- 
ply separated. The term g{p^) vanishes only for two 
cases: the trivial one which is free motion without dis- 
sipation/pumping, or the case of supercritical pump- 
ing where = Pq, eq. ( ^ ) for each particle. Then, 
the mean momentum becomes an invariant of motion, 
P{t) = Pq — const, and the center of mass moves ac- 
cording to R{t) = R{0) + Po{t). This behavior may also 
critically depend on the initial conditions of the particles, 
Pj(0), and shall be investigated in more detail now. 

In 1^ an approximation for the mean velocity P{t) of 
the swarm in one dimension is discussed which shows the 
existence of two different asymptotic solutions dependent 
on the noise intensity D and the initial momentum p^{0) 
of the particles. Below a critical noise intensity Dc, the 
initial condition p^ (0) > Pq leads to a swarm the center of 
which travels with a constant non-trivial mean velocity, 
while for the initial condition pf{Q) < Pq the center of 
the swarm is at rest. 

We can confirm these findings by means of two- 
dimensional computer simulations presented in Fig. ^ and 
Fig. ^, which show the mean squared displacement, the 
average squared velocity of the swarm and the squared 
mean velocity of the center of mass for the two different 
initial conditions. 

For (a) pf(0) > pg we find a continuous increase of 
AR^{t) (Fig. Ha), while the velocity of the center of mass 

reaches a constant value: P'^{t) = (^N^^ J^iPii^)^ ^ 

Pq, known as the stationary velocity of the force-free case, 
eq. (^). The average squared velocity reaches a constant 
non-trivial value, too, which depends on the noise inten- 
sity and the initial conditions, p|(0) > pg, i.e. on the 
energy initially supplied (cf Fig. H top) . 



10° 



10' 



< 

10" 



10" 



(a) 



(b) 



10" 



10' 



10^ 



10^ 



FIG. 5: Mean squared displacement AR^{t), eq. ^% 
of a swarm of 2.000 particles coupled to the center of 
mass. Initial conditions: rj(0) = {0.0,0.0}, (a) Pi{0) — 
{1.0, 1.0}, (b) p,(0) = {0.0,0.0}, for aU particles. Param- 
eters: a = l, D = 10"^ d2 = 10.0, So = 10.0, 70 = 20.0, 
c= 1.0. 




FIG. 6: Squared velocity of the center of mass, P^(t) = 
\ 2 

N^^ J2iPii^)) (solid hues) and averaged squared ve- 
locity p^(t) — N^^ Ylpfi^) (dashed lines) for the simula- 
tions shown in Fig. |[ (top) initial conditions (a), (bot- 
tom) initial conditions (b). 
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For (b) pf{Q) < Pq we find that the mean squared 
displacement after a transient time reaches a constant 
value, i.e. the center of mass comes to rest (Fig. || b), 
which corresponds to P'^{t) ^ in Fig. || (bottom). In 
this case, however, the averaged squared velocity of the 
swarm reaches the known stationary velocity, p'^{t) — 
1/N^pf{t) ^ Pq- Consequently, in this case the energy 
provided by the pumping goes into the motion relative to 
the center of mass while the motion of the center of mass 
is damped out (cf. also |^). Thus, in the following we 
want to investigate the relative motion of the particles in 
more detail. 

Using relative coordinates, {xi,yi} = ri — R, the dy- 
namics of each particle in the two-dimensional space is 
described by four coupled first-order differential equa- 
tions: 



Px: 



— Px 
— Px = ~l 



p. 

„2 



{pi)px^ - ax^ + {2Df'''i,{t) 



P, 



Py. -Py^ -.9 {pi)Py^ -ay^ + {2 0)^'^^^) (45) 

For P — Q, i.e. for the initial conditions pf < Pq and 
sufficiently long times, this dynamics is equivalent to the 
motion of free (or uncoupled) particles in a parabolic 
potential U{x,y) = a{x'^ + y'^)/2 with the origin {0,0}. 
Thus, within this approximation the system becomes a 
canonical-dissipative system again, in the strict sense 
used in Sect. ||. 

Fig. presents computer simulations of eq. ( |45| ) for the 
relative motion of the particle swarm in the parabolic po- 
tential. [43] (Note that in this case all particles started 
from the same position slightly outside the origin of 
the parabolic potential. This has been chosen in order 
to make the evolution of the different branches more 
visible.) As the snapshots of the spatial dispersion of 
the swarm show, we find after an inital stage the oc- 
curence of two branches of the swarm which results from 
a spontaneous symmetry break (cf. Fig. |^ top). These two 
branches will after a sufficiently long time move on two 
limit cycles (as already indicated in Fig. ^ bottom). One 
of these limit cycles refers to the left-handed, the other 
one to the right-handed direction of motion in the 2d- 
space. This finding also agrees with the the theoretical 
investigations of the deterministic case which showed 
the existence of a limit cycle with the amplitude 

-1/2 



ro = bo I 



(46) 



provided the relation so(i2 > 7oc is fulfilled. In the small 
noise limit the radius of the limit cycles shown in Fig. 
agrees with the value of rp. Further, Fig. || has shown 
that the averaged squared velocity p'^{t) of the swarm 
indeed approaches the theoretical value of eq. (|3^). 

The existence of two opposite rotational directions of 
the swarm can be also clearly seen from the distribution 
of the angular momenta Li of the particles. Fig. H shows 
the existence of a bimodal distribution for p{Lj. (The 
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FIG. 7: Snapshots (relative coordinates) of a swarm 
of 10.000 particles moving according to eq. (^) with 
P = 0. (top) t = 15, (bottom) t = 99. Initial conditions: 
{xi,yi} = {0.5,0.0}, {pxi,Pyi} = {0.0,0.0} for all parti- 
cles. Parameters: a = 1, D — 10^^, sq = 10.0; c = 1.0; 
70 = 20, d2 = 10. 



observant reader may notice that each of these peaks ac- 
tually consists of two subpeaks resulting from the initial 
conditions, which are still not forgotten a.t t — 99). Each 
of the main peaks is centered around the theoretical value 

m^Lo^roPa (47) 

where rg is given by eq. (Effl and po is given by eq. (^ 



The emergence of the two limit cycles means that the 
dispersion of the swarm is kept within certain spatial 
boundaries. This occurs after a transient time used to es- 
tablish the correlation between the individual particles. 
In the same manner as the motion of the particles be- 
comes correlated, the motion of the center of mass is 
slowed down until it comes to rest, as already shown in 
Fig. I 

This however is not the case if the initial conditions 
Pi{0) > Pq are chosen. Then, the energy provided by the 
pumping does not go completely into the relative motion 
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FIG. 8: Angular momentum distribution p{L) for a 
swarm of 10.000 particles at i = 99. The figure refers 
to the spatial snapshot of the swarm shown in Fig. |^ 
(bottom). 



of the particles and the establishment of the limit cy- 
cles as discussed above. Instead, the center of mass keeps 
moving as shown in Fig. ^ while the swarm itself does 
not establish an internal order. Fig. || displays a snapshot 
of the relative positions of the particles in this case (note 
the different scales of the axes compared to Fig. ^. 
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FIG. 9: Snapshot (relative coordinates) of a swarm of 
10.000 particles moving according to eq. (MS) at t = 
99. Initial conditions: {xi,yi} = {0.5,0.0}, {pxi,Pyi} — 
{1.0, 1.0} for all particles. Parameters see Fig. ^. 



B Coupling via Mean Momentum and Mean 
Angular Momentum 

In the following we want to discuss two other ways of 
global coupling of the swarm which fit into the general 



framework of canonical-dissipative systems outlined in 
Sect.||. There, we have introduced a dissipative potential 
G(/o, /i, /2, •■•) which depends on the different invariants 
of motion, J^. So far, we have only considered Iq — H, eq. 
( pT| ) in the swarm model. If we additionally include the 
mean momentum Ii = P, eq. (^4|) as the first invariant 
of motion, the dissipative potential may read: 



N 



G(/o,/i) = ^Go(p2) + Gi(P) 



Gi{P) 




NPi 



(48) 



(49) 



Here, Go(p^) is given by eq. (|36|). The stationary solution 
of the probability distribution, (p) is again given by eq. 
(p0|). In the absence of an external potential U{r), the 
Langevin equation that corresponds to the dissipation 
potential of eq. (E^) reads now: 



'9{p1)p,-Cp 
-{2DY'^Ut) 



(50) 



The term Gi (P) is choosen that way that it may drive 
the system towards the prescribed momentum Pi , where 
the relaxation time is proportional to Cp^. If we would 
have a vanishing dissipation function, i.e. 5 = for pf = 
Pq, it follows from eq. (50) for the mean momentum: 



N 



N 



C-pt 



(51) 



The existence of two terms Go and Gi however could 
lead to competing infiuences of the resulting forces, and 
a more complex dynamics of the swarm results. As before, 
this may also depend on the initial conditions, i.e. > 

Po or Pl< pI- 

Fig. |l^ shows the squared velocity of the center of 
mass, 

p2(i) and the averaged squared velocity of the 
swarm p^{t) for P^ < Pq, We find an intermediate stage, 
where both velocities are equal, before the global coupling 
drives the mean momentum P towards the prescribed 
value Pi, i.e. P^{t) ^ (P^^ + P^y). On the other hand, 
p^it) — > Pq, as we have found before for the force- free 
case and for the linearly coupled case for similar initial 
conditions. The noticeable decrease of P^ after the initial 
time lag can be best understood by looking at the spatial 
snapshots of the swarm provided in Fig. 11. For t = 10, 
we find a rather compact swarm where all particles move 
into the same (prescribed) direction. For t = 50, the cor- 
relations between the particles have already become ef- 
fective, which means the swarm begins to establish a cir- 
cular front, which however does not become a full circle 
[44]. Eventually, we find again that the energy provided 
by the pumping goes into the motion of the particles rel- 
ative to the center of mass, while the motion of the center 
of mass itself is driven by the prescribed momentum. 
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FIG. 10: Squared velocity of the center of mass, P^{t) — 
(^N-^J2tPi(t)) (solid lines) and averaged squared ve- 
locity p^{t) — N^^ J2Pi(^) (dashed lines) for the simula- 
tions shown in Fig. 11. 

For the initial condition > the situation is differ- 
ent again, as Fig. |l2| shows. Apparently, both curves are 
the same for a rather small noise intensity, i.e. P^{t) = 
p^(t) are both equal, but different fromp§, eq. and 
the prescribed momentum Pf. This can be only realized 
if all particles move in parallel into the same direction. 
Thus, a snapshot of the swarm would much look like the 
top part of Fig. 11. 

Eventually, we may also use the second invariant of mo- 
tion, I2 = L, eq. (Eq) for a global coupling of the swarm. 
In this case, the dissipative potential may be defined as 
follows: 



N 



G{Io,l2) = Y.Go{pl) + G2{L) 



(52) 



Cr 



G2{L) = ^r£r,xp,-NL, 



(53) 

G'o(p^) is again given by eq. ( |36|) . The term 6*2 (i) shall 
drive the system to a prescribed angular momentum Li 
with a relaxation time proportional to C^^ . 

Li can be used to break the symmetry of the swarm to- 
wards a prescribed rotational direction. In Sect. IV A we 



have observed the spontaneous occurence of lefthand and 
righthand rotations of a swarm of linearly coupled parti- 
cles. Without an additional coupling, both rotational di- 
rections are equally probable in the stationary limit. Con- 
sidering both the parabohc potential U{r,R), eq. ( pl| ) 
and the dissipative potential eq. (|5^), the corresponding 
Langevin equation may read now: 



N 



Cl r, X [S^ Vi X p^- N 
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FIG. 11: Snapshots of a swarm of 2.000 particles moving 
according to eq. (|o|). (top) t = 10, (bottom) t — 50. 
Initial conditions: {r„,ryi} = {0.0,0.0}, {pxi,Pyi} = 
{0.0,0.0} for all particles. Parameters: {Pix,Piy} = 
{0.344,0.344}, Cp = 10'^ D = 10-^ sq = 10.0; c = 1.0; 
70 = 20, d2 = 10. 

+ i2D)'/'Ut) (54) 

The computer simulations shown in Fig. |l^ clearly dis- 
play an unimodal distribution of the angular momenta Li 
of the particles, which can be compared to Fig. || with- 
out coupling to the angular momentum. Consequently, 
we find in the long time limit only one limit cycle corre- 
sponding to the movement of the swarm into the same 
rotational direction. The radius tq of the limit cycle is 
again given by eq. (|46|). 

We would like to add that also in this case the dy- 
namics depends on the initial condition, Li. For simplic- 
ity, we have assumed here |Xi| = Lq = ropo, eq. ( |47| ) 
which is also reached by the mean angular momentum L 
in the course of time (cf . Fig. |l^) . For initial conditions 
|ii| ^ Lq, there is of course no need for the rotation of 
all particles into the same direction. Hence we observe 
both left- and righthanded rotations of the particles with 
different shares, so that the mean angular momentum 
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FIG. 12: Squared velocity of the center of mass, 
P^{t) = (^N~^J2iPiit)^ (solid lines) and averaged 

squared velocity = N^^^p'j{t) (dashed Hnes) for 

a swarm of 2.000 particles moving according to eq. (50). 
{Pix, Ply} — {1-0, 1.0}, for the other parameters and ini- 
tial conditions see Fig. 11. 

is still i( — > ii. This results in a broader distribution 
of the angular momenta of the particles instead of the 
clear unimodal distribution shown in Fig. |l^. For initial 
conditions \L\ \ ^ Lg on the other hand, the stable rota- 
tion of the swarm breaks down after some time, since the 
driving force L Li tends to destabilize the attractor 
L ^ Lq. This effect will be investigated in a forthcom- 
ing paper, together with some combined effects of the 
different global couplings. 



V. DISCUSSION 

Eventually, we can also combine the different global 
couplings discussed above by defining the dissipation po- 
tential as: 

G{Io,Ii,l2) = G(p^P,L) 

^ Goip^) + GiiP) + G2iL) (55) 

Go(p') isgiven by eq. (||), Gi(P) by eq. ® and G2(i) 
by eq. (|53| ) . Considering further an additional - external 
or interaction - potential, the corresponding Langevin 
equation can be written in the more general form: 



Pi 



d 



Pi = ~g{p^ )Pt -ir-io^oU (r, R)] 



Or, 

d 

dp. 



ai (P-Pi)Va2 {L^LiY 



(56) 



FIG. 13: Angular momentum distribution p{L) for a 
swarm of 2.000 particles at t = 99. For comparison with 
Fig. ^ the parameters, initial conditions and snapshot 
times are the same as in Fig. 0. Additional coupling: 
Cl - 0.05,, {Li J - {0.4}. 



The mean momentum P and mean angular momentum 
L are given by eqs. (p^), (p5[); whereas the constant vec- 
tors Pi and Li are used to break the spatial or rota- 
tional symmetry of the motion toward a preferred di- 
rection. The different constants ai may decide whether 
the respective influence of the conservative or dissipative 
potential is effective or not, they further determine the 
time scale when the global coupling becomes effective. 
The term g{p^), eq. ( p9| ) on the other hand considers the 
energetic conditions for the active motion of the swarm, 
i.e. it determines whether the particle of the swarm are 
able to "take off" at all. 

The combination of the different types of coupling may 
lead to a rather complex swarm dynamics, as already in- 
dicated in the examples discussed in this paper. In par- 
ticular, we note that the different terms may have com- 
peting influences on the swarm, which then would lead to 
a "frustrated" dynamics with many possible attractors. 

In this paper, we have basically restricted the investi- 
gation of the swarm dynamics to global couplings which 
fit into (or can be reduced within some approximations 
to) the general outline of canonical-dissipative systems. 
Finally, we want to add some comments on that. On one 
hand, it is possible to extend this kind of approach also 
to other invariants of motion, this way e.g. covering pre- 
vious investigations of swarms coupled via the mean ori- 
entation of the particles [Q, ||. On the other hand, we 
want to underline that canonical-dissipative systems are 
a theoretical class of models, where both conservative and 
dissipative elements of the dynamics are determined by 
invariants of the mechanical motion. Thus, from this per- 
spective, a more realistic swarm dynamics may be also 
based on less restrictive assumptions. 

The advantage to use canonical-dissipative systems as 
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a framework for swarm dynamics is given by the fact 
that in many cases the rather complex dynamics can 
be mapped to an analytically tractable model. With the 
Hamiltonian theory of many-particle systems as a start- 
ing point, we are able to extend known solutions for con- 
servative systems to non-conservative systems. This al- 
lows us to to construct a canonical-dissipative system, 
the solutions of which converges to the solution of the 
conservative system with given energy. That means, for 
given initial conditions, we could predict the asymptotic 
solution of the canonical-dissipative dynamics by means 
of the solutions of the Hamiltonian equations on the re- 
spective energy surface. 

In addition to these theoretical considerations which 
have a value of their own, we want to note that the frame- 
work of canonical-dissipative systems still covers impor- 
tant features of real (biological) systems, such as energy 
take-up and dissipation. The general description outlined 
in this paper allows us to gradually add more and more 
complexity to the swarm model, this way bridging the 
gap between a known (physical) dynamics to a more com- 
plex (biological) dynamics Some hints for this 
shall be given at the end. 

On the level of the "individual" particles, the key dy- 
namics of the model is given by a modified Langevin 
equation that, in addition to stochastic influences, also 
considers other forces on the particle, resulting e.g. from 
external potentials, interactions, use of stored energy [^ , 
or influences of a self-consistent field ||^, |j that already 
exceed the framework of canonical-dissipative systems. 

The consideration of an external potential also allows 
to model the spatial environment of the swarm, for in- 
stance to consider obstacles |30 . Additionally, we can 



also consider that the pumping of energy for the parti- 
cles is restricted to certain spatial domains which model 
food sources. In this case the dissipation function g{H) 
also becomes a spatial function. Such an extension has 
been already discussed in |^0|, ^ and can be also imple- 
mented in the swarm model discussed here. Eventually, 
we note that the genuin particle-based approach to col- 
lective phenomena used in this paper is not restricted to 
biological systems, but also applicable for describing and 
simulating complex interactive systems in a wide range of 
applications, even in economics and social systems [M. 
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